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ABSTRACT 




£ 



The basic filter-observer equations of Kalman for optimal and 
suboptimal filters are studied using the concepts of Lyapunov functions 
and stability theory. The Second Method of Lyapunov is used to form a 
basis for comparison of the convergence rates of such filters. Lyapunov 
functions are also used to derive constraining relations for the 
elements of the filter gain matrix leading to design criteria for sub- 
optimal filters. A derivation of the optimal filter gain based upon the 
Lyapunov function of a random variable is given. This derivation shows 
that the optimal filter converges most rapidly. A design of a sub- 
optimal filter for one class of signal models is given based solely upon 
stability constraints. 
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SYMBOLS AND ABBREVIATIONS 



X 

x 

t 

det(A) 

tr(A) 

X(A) 

INI 

l“l 

N(x e ,6) 

E[x] 

£ 

inf [ * ] 



an n square matrix unless otherwise specified 

an n vector. Scalars will be specified unless obvious 
from the notation 

independent variable time and is always a scalar 
(similarily for n in discrete time) 

determinant of matrix A 

trace of matrix A 

eigenvalue of matrix A 

norm of the argument (see Appendix A) 

magnitude of the complex scalar a 

neighborhood about the point x g of radius 6 . 

Taken as the usual sphere definition all x such that 
I |x - x e | | i 6 

expected value of the random vector x 
is a member of, belongs to 
greatest lower bound 
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I. INTRODUCTION 



For many problems in control it is desirable to have some measure 
of all states in a given signal process. By their nature, many 

processes only allow us to physically measure or observe some of the 
states characterizing the process. Moreover, there is some amount of 
uncertainty in the measured states which may be characterized by 
additive measurement noise. Thus, estimation theory, as a device for 
estimating the states of a system in the presence of noisy and incom- 
plete observations, has emerged as an important part of modern control 
theory. 

In the early 1 960 1 s Kalman [K3,K4] advanced the well known Wiener 
filter theory [ref. D2,P1] by showing how such linear least square 
filters, which Wiener synthesized using classical frequency domain 
theory, may be realized in the time domain. This development coin- 
cided with the development of the digital computer and its use as an 
active component in process control systems. 

The basic filter equation of Kalman is really a recursive weighted 
average technique of the following form: 

X = X + G (X-Z) (1) 

where X is the corrected estimate, X is the projected estimate 
following known signal dynamics, G is the filter gain weighting 
factor, and Z is the signal observation including noise. G is 
generally time varying and its selection constitutes the basic design 
and nature of the filter. The difference equation of (1) becomes a 
differential equation in continuous time. 
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The corrected estimate, X , is usually an n vector of the states 
of the signal process, and the observation, Z , is an m vector with 
m <_ n. 

The selection of the filter gain has two objectives: 

1. To provide an estimate of unobserved states 
based upon data from the observed states. 

2. To provide for some weighting of the observations 
to allow compensation for the measurement noise. 

Objective 1 is the basic observer problem which has emerged in recent 

years (such as the Luenberger observer [LI]) and does not consider 

any effects of noisy measurements. Objective 2, as stated, involves 

estimation theory and the stochastic properties of the observation 

noise as well as the random input signals which drive the system. 

Current literature contains many examples of the application of 
Kalman filters to various problems. The computational capability 
required to implement such filters is general ly “large, and often a 
much simpler implementation works nearly as well. The literature also 
indicates that the method generally employed to determine the parame- 
ters of these filter implementations involves large scale simulation 
trial s. 

This thesis studies the properties of the basic filter observer 
equation in both continuous and discrete time using the Second Method 
of Lyapunov. Lyapunov theory [ref.. K1,H1] was developed for investi- 
gating the stability properties of systems. Its application to the 
filter-observer problem has appeared in recent literature (Deyst and 
Price [D1 ] ) . However, this approach has not been fully exploited in 
the literature, and its possibilities are investigated in this thesis. 
This study leads to engineering insight into the operation of suboptimal 
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filters and provides some basis for comparison of the performance of 
these filters when applied to the same problem. Some specific results 
obtained include constraining relations on the filter gain matrix 
which insure: 

1. that the filter is asymptotically stable. 

2. that the filter will converge to steady state 
at a rate faster than a given bound. 

Application of the filter observer equations as discussed here is 
limited to linear systems and it is assumed that the signal dynamics 
are known exactly. A review of stability concepts for linear time- 
varying systems and the Second Method of Lyapunov is included. 

Chapter II reviews basic definitions and stability theorems for 
linear' systems, based upon the concepts of state-space theory, with 
emphasis on forced systems so that these results may be applied readily 
to the filter-observer problem. 

Chapter III reviews the Second Method of Lyapunov. This chapter 
includes some theorems on the determination of Lyapunov functions for 
certain classes of systems as well as some theorems on transient esti- 
mation using Lyapunov theory which have not appeared in the literature 
previously. 

Chapter IV applies stability theory to the filter-observer homo- 
geneous dynamics and contains many original results including the 
formulation of a time invarient Lyapunov function for the filter. The 
same Lyapunov function may be used for a large class of filters 
applied to the same problem. Particular filters have different time 
derivatives of the Lyapunov function. This leads to a method of 
comparing convergence rates for various filters applied to the same 
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problem. By applying the Hadamard-Gerschgorin theorem to provide 
sufficient conditions for positive definiteness of a symmetric matrix, 
certain constraining relations are derived for the elements of the 
filter gain matrix which lead to the design of stable filters. 

Chapter V considers the forced filter dynamics. New results 
include the introduction of the concept of a Lyapunov function for a 
random variable. A derivation, using these Lyapunov functions, of 
the optimal filter is given showing that such filters also converge 
most rapidly to the minimum covariance of error. 

Chapter VI applies the results of Chapter IV to the design of a 

J_ L. 

class of suboptimal filters for the general n u order signal model 
in phase variable form. 
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II. STABILITY OF LINEAR SYSTEMS 



A. INTRODUCTION 



A review of stability definitions and theorems are set forth in the 
following., Many current references are available for this material 
(HI, K1 , K2, 01, SI, T1 and others). Although here we are mainly 
interested in linear systems, many of the theorems, as noted, apply 
to nonlinear systems. The reader is referred to Appendix A for a 
summary of norms and norm properties. 



Definition 2.1 Free system: Any system with no input 

forcing functions. Such a system is 
described as x = f(x,t) 



De finition 2.2 Autonomous System: A free system which 

is also time invariant, (i.e., x = f(x)) 



Definition 2.3 



Definition 2.4 



Definition 2.5 



Equilibrium state(x e ): Any state for 

which x = 0. For a free system 
f(x e ,t) = 0 

Boundedness: An equilibrium state x g 

is said to be bounded if, and only if, 
there exists a neighborhood N(x ,6) 
such that, if the initial state x(t Q ) 
lies within N, then l|x(t)-x e || <_ m < °° 
for all t > t 

o 

Stability (in the sense of Lyapunov): 

An equilibrium state x g is said to be 
stable, if, and only if, to any neigh- 
borhood N(x e ,e) there" corresponds a 
neighborhood N(x e ,<$) such that, if 
x(t Q ) lies in N(x g ,6) then x(t) lies 

in N(x ,e) for all t > t . 
e o 
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The difference between the last two definitions is important for 
nonlinear systems although the definitions are equivalent for linear 
systems. Note that Definition 2.5 requires that the state remain 
within a' preassigned neighborhood (which may be arbitrarily small) for 
all t > t whereas the definition of boundedness requires only that the 
distance (in state space) from x g to the trajectory x(t) remain bounded. 
Thus, boundedness allows for such things as limit cycles. 

Definition 2.6 Asymptotic Stability: An equilibrium state 

is said to be asymptotically stable if 

(1) It is stable. 

(2) | |x(t)-x e | | approaches zero as t 
approaches infinity for any initial con- 
dition lying in the neighborhood N(x e ,6) 
for which it is stable. 

The above definitions are for strictly local conditions about the 
equilibrium state. If the neighborhood N(x e ,<$) can be the entire state 
space, then the system is said to be globally stable or asymptotically 
stable in the large (ASIL) about the equilibrium state x J 

Uniform stability is also a useful concept. It means that the 
neighborhood in which the initial state must lie is independent of the 
initial time. That is, <$ is not a function of t . 

Definitions 2.1 through 2.6 above also apply to discrete time 
systems with the use of (n,n Q ) in place of ( t , t Q ) respectively. In 
the sequel, theorems and concepts are developed for both continuous 
and discrete time systems. The discrete notation is as follows: 



^Some authors distinguish between global and stability in the 
large (HI) but for linear systems, they are identical. 
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( 2 . 1 ) 



The System : 

x(n+l) = A D (n) x(n) + B Q (n) u(n) 
The Observation : 

y(n) = H D (n) x(n) 

The Fundamental Matrix: 



*(n,n Q ) = A D (n-l) A Q (n-2) . . . Ap(n Q ) 



The Solution: 



n-1 



x(n) = $(n,n ) x(n ) + e $(n,k+l) B(k) u(k) 

0 0 k=n 



( 2 . 2 ) 

(2.3) 

(2.4) 



For time invariant systems we have $(n) = A^ n . 

Also it is interesting to note the relationship between continuous 
sampled systems and their discrete equivalent, for the time invariant 

case. 

T 

A d = $(T) B D = / $(T-t) B dx (2.5) 



B. FREE SYSTEMS 

It is well known that for linear time invariant systems, if the 
eigenvalues of the system matrix A all have negative real parts, then 
the system is asymptotically stable. Similarly for discrete systems 
the system is asymptotically stable to the origin if, and only if, the 
magnitude of the eigenvalues of Ap are less than unity. Such condi- 
tions are not easily established for time varying systems. However, 
the following theorems apply for the system x(t) = A(t) x(t), with 
x g = 0 and the equivalent discrete case [ref. SI]. 

T heorem 2-1 The origin of a linear system is bounded 

if, and only if, the norm of the fundamental 
matrix $(t,t Q ) (or 3>(n,n Q )) is bounded. 
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To prove sufficiency note that 

||x(t)|| = | |*(t,t 0 ) x(t 0 ) II 1 ||*(t.t 0 )|| ||x(t 0 )|| 

If ||4(t,t 0 )|| <( 

Then | |x(t)| | £ K| |x(t Q )| | 

By definition 2.4, with 6 = | | x( t Q ) | | it follows that the origin 

( x e =0) is bounded. Necessity is shown by considering 
| | x( t) | | = | |*(t,t Q ) x(t Q )| | 

and noting that if the norm of $(t,t Q ) is not bounded then some 

element of $(t,t Q ) approaches infinity with time and hence ||x(t)|| 

also approaches infinity and is unbounded. 

Theorem 2.2 Boundedness and stability of the 

origin are equivalent in linear systems. 

To prove this, we must show that for any e > 0 it is possible to find 

a 6 such that if x(t Q ) is in N(0 ,6 ) then x(t) is in N(0,e) for all 

t > t Q . This reduces to showing that if | | x( t Q ) | | < 6 then ||x(t)|| < e. 

By definition 2.4, if the origin is bounded we may take 6 = . 

Then by theorem 2.1 

I |x(t)| I < K||x(t 0 )|| <-^f- = e 

which shows that boundedness implies stability. It is obvious that 

stability implies boundedness hence the two are equivalent. 

Theorem 2.3 The origin of a linear system is 

asymptotically stable if, and only if, 
the norm of $(t,t Q ) is bounded for all 
t < t and 4>(t,t o ) approaches zero as 
t approaches infinity. 
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If ||$(t,t || <_ K then the origin is stable (Theorems 2.1 and 2.2). We 
now have to show that ||x(t)|| approaches zero as t approaches infinity 
Now | | x ( t ) | | = ||$(t,t Q ) x(t Q )|| and if $(t,t Q ) approaches zero then 

| | x ( t) | | approaches zero. 

Moreover, if ||x(t)|| approaches zero for all 1 1 x ( t Q ) 1 1 in N(0,6) 
then o(t,t 0 ) approaches zero, since x(t Q ) is arbitrary in N(0,6). 

It should be noted that if a linear system is asymptotically stable 
then it is also ASIL, since the conditions for stability depend only 
upon the fundamental matrix and are independent of the initial condi- 
tions of the states. 

The foregoing theorems also apply to discrete systems with the 
change in notation from $(t,t Q ) to $(n,n Q ). 



C. FORCED SYSTEMS 



Now we are interested in studying forced systems of the form 

x(t) = A(t) x(t) + B(t) u(t) (2.6) 

We are particularly interested in the effect of u(t) on the system, 
assuming it is at rest at t=t . Then the solution of (2.6) is 



x(t) = / $(x,t ) B(t) u(t) dt 

*0 



(2.7) 



Stability of forced systems is considered with respect to a given 
set of inputs U = { u ( t ) > . Generally, this set is taken to be the 
set of bounded inputs. 

Definition 2.7 Stability of a forced system. 

A forced system is stable with 

respect to U if, and only if, 

the states are bounded (||x(t)|| 

<_ K) for all u(t) in U and all 

t > t . 

— o 
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Theorem 2.4 A forced system is stable with 

respect to the set of bounded inputs 

1f t 

l I | <*» ( t , t ) B ( x ) | | dt <_ K, 
x o 

To prove sufficiency note that 

t t 

I I X ( t) 1 1 1 / I |$(T,t )B(t)u(t)| I dt <_ / | |*(t,t )B(t)| | I I u (t ) j I dt 

t„ ° t 

O o 

If I I u(x) 1 1 £ K 2 then | |x(t)| | < K^ 2 . 

Hence the state is bounded. 

A consolidating theorem appears in Timothy and Bona [Tl] and 

Kalman [Kl] with proof given in [K1 ] . 

Theorem 2.5 For the system described by equation 

(2.6) if 

(i) I |A(t)| | <. K 1 for all t 

(ii) 0 < K 2 < | | B ( t ) v | | < K 3 for all ||v||=l 

where v is an arbitrary vector, 
then the following statements are 
equivalent. 

1. Any uniformly bounded input 
1 1 u (t) 1 1 <_ K^ results in a 
bounded output ||x(t)|| < Kg < “ 

2. For all t > t 

o 

t 

/ | |<t>(t,T) | |dtl K g < » 

t o 

3. The free system is uniformly 
asymptotically stable. 
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III. THE SECOND METHOD OF LYAPUNOV 



A. INTRODUCTION 

The foregoing indicates that stability for linear systems is 
completely determined by the boundedness properties of the fundamental 
matrix and the inputs to the system. This means in essence that the 
complete solution to the system must be obtained before stability can 
be ascertained. Lyapunov proposed a method by which stability of a 
system may be determined without finding the solution to the system 
equations. His method is based upon a generalized energy concept. 

From physics it is known that for a nonconservative system, if the 
total energy is always decreasing, then the system is stable. The 
total energy of such a system is a scalar function of the state of the 
system. Lyapunov has shown that the very existence of scalar functions 
of the state of a system which obey certain properties is sufficient 
to conclude that the system is stable. These functions are called 
Lyapunov functions. It is interesting to note that Lyapunov functions 
for a stable system are not unique and that a sufficient condition for 
stability is the existence of any such function. 

Definition 3.1 Lyapunov function. Any function, V(x), 

having the following properties is 
called a Lyapunov function. 

(a) V(x) is a continuous scalar function 
of the systems state vector and has 
continuous first partial derivatives. 

(b) V(x) is positive definite. 

(c) V(x) = [grad V(x)]^ x is negative 
defi nite. 
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Definition 3.2 A scalar function, V(x), of a vector argu- 
ment is positive definite if, and only if, 

(a) V(0) = 0 

(b) V(x) > 0 whenever x f 0 

Positive semi -definite means that equality 
is also allowed in part (b). 

For negative definite functions the 
inequality is reversed. 

In two dimensional state space the Lyapunov function is a cupped 
shaped surface setting on the origin. Its continuity and definiteness 
mean that the level curves of V projected on the state plane are 
closed' about the origin and the curve = V (x ) is inclosed within 
k 2 = V(x) if k^ < k^ . 

Under the foregoing conditions we can consider that the Lyapunov 

2 

function gives a measure of distance from the origin in state space 
as the system follows the trajectory x(t). 

B. THE BASIC STABILITY THEOREM 

Before stating the theorem, consider the effect of explicit time 
variation upon the properties of the Lyapunov function. Suppose we 
have a function V(x,t) and its time derivative 

»(x,t) = -ff- + [grad V] T x (3.1) 

Now, V(x,t) will be positive for all t if it is always greater 
than a positive function e(||x||). However, this is not sufficient 
to guarantee asymptotic stability even if < 0 for all t and x f 0. 



2 

Distance used in this sense is not the same as the usual 
Euclidian distance. Distance is implied by the ordered nature of the 
level curves of V as discussed. 
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Since the grad V term in equation (3.1) indicates the motion of the 
system state (i.e., is multiplied by x), it is possible that 
be sufficiently negative so that < 0 for all t, while the state 
moves outside of any bounding region (N(x g »e) ) . Hence, the state may 
never go to the equilibrium state x g . Hsu and Meyer [HI] have a good 
discussion and Kalman and Bertram [Kl] provide a rigorous proof of 
the following theorem. 

Theorem 3.1 Consider the continuous time, free 

dynamics system x = f(x,t) where f(0,t) = 0 
for all t. Suppose there exists a scalar 
function V(x,t) with continuous first 
partial derivatives such that V(0,t) = 0 and 

(i) there exists a continuous, non- 
decreasing scalar function a such that 
a(0) = 0 and V(x,t) >_ ct( | | x 1 1 ) for all t 
and x f 0. 

(ii) there exists a continuous positive 
scalar function y such that y(0) = 0 and 
¥(x,t) < - Y (||x||) 

(iii) V(x,t) < e( | |x 1 1 ) 

where e is a continuous nondecreasing 
function such that g(0) = 0 for all t. 

(iv) a[ 1 1 x 1 1 ) -* * with | |x| | * 

Then the system is asymptotically stable in the large to the 
equilibrium state x g = 0 and V(x,t) is a Lyapunov function for the 
system. 

Requirement (i) insures that the Lyapunov function is always 
positive whereas (iii) insures that the function does not stay infinitely 
large throughout the state space. In particular it is required that the 
V function doeskin fact, go to zero uniformly with time at the origin. 
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The y function in (ii) assures that V is always negative. This 
together with (i) and (iii) assure that the Lyapunov function following 
the system trajectory is decreasing and must end up at the origin. 

To insure that the conditions exist for all of state space, require- 
ment (iv) is imposed and the system is ASIL. It should be noted that 
a positive definite quadratic form for V automatically satisfies the 
requirements of (i), (iii), and (iv). 

For discrete systems we define a Lyapunov function as V(x,n) 
and the rate of increase of V along the system trajectory as 

AV(x,n) = V(x(n+1 ) ,n+l ) - V(x(n),n) (3.2). 

With these changes, Theorem 3.1 becomes the corresponding theorem 
for discrete systems. 

C. THE LINEAR AUTONOMOUS SYSTEM 

For linear autonomous systems, it is a necessary and sufficient 

condition for ASIL that a quadratic form Lyapunov function exist. 

This is embodied in a theorem originally proved by Lyapunov. 

Theorem 3.2 - The origin of a linear autonomous 

system, is asymptotically stable if, 
and only if, for any symmetric positive- 
definite matrix, C, there exists a 
unique positive-definite matrix Q which 
satisfies the matrix equation 

A T Q + QA = -C (3.3) 

A full proof of this theorem appears in Hahn [H2]. Essential parts 

also are shown by Bellman [B1 ] . The major results of interest are 

delineated here for convenience. 

1. The Lyapunov function is V = x T Qx 
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2. By direct computation , 

= x T Qx + x T Qx = x T A T Qx + x Vx 

t (3.4) 

= -x Cx 



3. 



4. 



Given any positive definite matrix C then Q is uniquely 
determined by the matrix equation (3.3) if no two eigen- 
values of A sum to zero or no eigenvalue is zero, (see [B1 ] . ) 
If A is stable then the unique Q is given by 



Q = 



” _A T t r At 
/ e C e dt 

0 



(3.5) 



Using 4 and a theorem 0 stated in Browne [B2] it is easy to show 



that if C is positive definite, then 0 must be also. Since C = B B 

3 

where B is a real non-singular matrix we can write 



Q = / e 

0 



ATt B T Be At dt = / (Be At ) T (Be At ) dt 



(3.6) 



At At 

But B and e are non-singular. Hence Be is non-singular and Q is 

3 

positive definite by reapplication of Browne's theorem . 

Generally the use of Theorem 3.2 proceeds as follows: 

(i) Choose any matrix C which is positive definite and solve the 
n(n+l)/2 algebraic equations (3.3) for the elements of Q. 

(ii) Test Q for definiteness. If it is positive definite the system 
is asymptotically stable. If' it is not positive definite the system 
is not asymptotically stable. This follows from the necessary and 
sufficient conditions stated in Theorem 3.2. 



3 T 

From Browne; If C is any real non-singular matrix and A = C C 

then A is positive definite. The converse is true, any positive 

definite matrix A can be expressed as such a product. Similarily, if 

C is n square of rank r<n then A is positive semi -definite of rank r. 
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The theorem corresponding to Theorem 3.2 for discrete time systems 
is as follows [ref. K1 ,01 ] . 

Theorem 3. 2D The discrete time system described by 

x(n+l) = A d x(n) 

is asymptotically stable if, and only if, 
for any given positive definite matrix C 
there exists a unique positive definite 
matrix Q satisfying 

A d T QA d - Q = -C (3.7) 



where V = x^Qx is a Lyapunov function 
with A V = -x T Cx. 

An algorithm for numerical solution of equation (3.7) for Q given 
any positive definite matrix C is derived in Appendix B. This is a 
new result and makes the application of Theorem 3. 2D much easier than 
in the past. 

At this point it is important to note that we cannot in general 
choose any positive definite matrix for Q and then perform the multipl 
cation and addition indicated in (3.3) to obtain a C matrix which is 
positive definite, even if A is known to be stable. The following 
example illustrates this. 

Example 1 . For A = 




it follows that 
C = -(A T + A) 



0 1 

1 4 



The system is stable, but C is indefinite. Therefore, x T x is not a 
Lyapunov function for the system described by the given transition 
matrix. 
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The following theorem has not generally appeared in the literature 
and is formulated here as a new result. 

Theorem 3.3 If, for a stable system, its transition 

matrix A is real and commutes with its 
transpose (i.e., AA^ = A^A) then 

V = | |x| | 2 (i.e., Q = I) 

is a Lyapunov function for the system 
described by A. 

This theorem gives a sufficient condition for the square of the 
norm of the states to be a Lyapunov function. To show this we take 
| |x| | = (x T x) 1/2 hence V = x^x and show that A*'" + A must be 

negative definite. But A^ + A is a real symmetric matrix and must 
have real negative eigenvalues [B2] if it is to be negative definite. 

But this is to say that the linear system described by the transition 
matrix A T + A is stable. Its fundamental matrix is then 
$ = . Now if A T A = AA T then we can write $ = e^ t e^. 

By Theorems 2.1 and 2.2 of Section I IB, if ||$|| <_ k then the system 

is stable. But since A is stable ||e | | <_ k.| so we write 

I I'll ik,k, = k • 

For discrete systems the following theorem applies. 

Theorem 3. 3D If A^ A^ = Aq A^ then a 

Lyapunov function for the discrete 
system x(n+l) = A^ x(n) is 

V = x T x with aV = x T (A^ Ap - I) x 

Suppose there exists some transformation of states x = Sy 
such that the transformed system matrix D - S _1 AS is normal 
(i.e., D T D = DD T ). Then we can apply Theorem 3.3 and a Lyapunov 
function for the transformed system is = y^y. Since stability 



25 



properties for linear systems are invariant under a similarity trans- 
formation, we may apply the inverse transformation to to find a 
suitable Lyapunov function in the x state space. The transformation 
must also be applied to V. This is summarized in the following theorem. 

Theorem 3.4 If a non-singular transformation 

exists for x = S _1 y such that 
D = SAS"^ is normal, then a 
Lyapunov function for the system 
x = Ax is: 

V = xVsx with 

V = x T S T (D T + D) Sx 

= x T (A T S T S + S T SA) x 

A similarity transformation may be applied to discrete systems as 
in the development of Theorem 3.4 for continuous systems. 

D. LINEAR FREE SYSTEM 

Here we look briefly at a system of the form x = A(t) x . 

According to Theorem 3.1 of this chapter, if we can find a Lyapunov 
function for the system it may be time varying. However, we may also 
be able to find a Lyapunov function which is not explicitly a function 
of time. Assuming that we can find a positive definite 0 such that 
V = x"*Qx is a Lyapunov function then 

v- x T [A^(t)Q + QA( t ) ] x = -x T C(t) x 
In order to apply Theorem 3.1, we must find a positive function 
y( | | x | | ) such that 

V(x,t) <_ -y(||x||) for all t. 

Apparently we must somehow remove the time variation of A(t). 
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This obviously depends upon the nature of A(t). A general pursuit of 
this problem is not very fruitful. However, for the special case 
considered next we can get some useful results. 



E. LINEAR SYSTEM WITH PERTURBED COEFFICIENTS 

The special case of the linear free system considered here is that 
the state transition matrix A(t) is of the form A(t) = A + G(t). 



Now we assume V(x) = x Qx and compute 
^(x.t) = -x T Cx + x T (G T (t)Q + 0G(t))x 



(3.9) 



where -C = A Q + QA . 

From the results in section 3C, we know that if the system x = Ax 
is stable then for any choice of a positive definite matrix C we can 
compute a positive definite matrix Q. Moreover, if G(t) tends to zero 
as t increases, then there exists some t such that for every t > t 
the first term in the expression for Hx,t) will dominate and is 
negative definite. Hence, the assumed quadratic form for V(x) is 
indeed a Lyapunov function. In this case we do not have uniform stabil- 
ity since the result depends upon t Q . 



Example 2, 



Let A(t) = 



-1 + 



-2 + 





-1 0 




1 


0 


Then A = 




G(t) = 


t 


i 

t_ 




0 -2 




0 



Applying Theorem 3.3 page 25 we take Q = I then C = 



2 0 
_0 4_ 

and C will dominate for all t >_ 1/2. Thus this system is ASIL for 
any t >_ 1/2. 
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For the discrete time system 

x (n + 1) = Ax(n) + G(n)x(n) 

we again assume a Lyapunov function V = x^Qx. Then by direct compu- 
tation 

A V = V(x(n+1)) - V(x(n)) = -x T (n)Cx(n) + x"^(n)QG(n)x(n) 

where -C = A^QA - Q . 

From Theorem 3.2, page 22, we know that for any positive definite 
C we can find a positive definite Q if A is stable. Similar comments 
apply if G(n) tends to zero as for the continuous case. 

F. FORCED LINEAR SYSTEMS 

The usual application of Lyapunov's Second Method is to the study 
of the dynamics of homogeneous linear systems or to the study of the 
dynamics describing the motion of a system about some nominal trajectory 
usually attributed to the forcing function. In this section we consider 
a different approach which leads to a novel interpretation of the 
Lyapunov function with inputs present. 

Consider the linear system 

x = Ax(t) + Bu(t) (3.10) 

Again assume a Lyapunov function V = x T Qx and compute its time 
derivative. 

. .T T . 

V = x Qx + x Qx 

= (u T B T + x T A T )Qx + x T Q(Ax + Bu) 

= x T (A T Q + QA)x + u T B T Qx + x T QBu 

V = -x T Cx + 2 x T QBu (3.11) 

From Section IIC we know that if the system x = Ax is stable we 
can find a positive definite Q for any given positive definite C such 
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that -C = A^Q + QA. Also, note that the bilinear forms u"Vqx and 
x T QBu are transposes of each other. Since they are scalars, they must 
be equal, and they may be combined as shown in the second term of 
equation (3.11). 

By definition 2.1 in Section 1 1 C page 13 the system described by 
(3.10) is stable if and only if ||x|| £ . Thus if we can show that 

0 given by (3.11) is negative for ||x|| > k.] we can be assured that 
the Lyapunov function will converge to a region determined by ||x|| <_ k y 
This result is stated in the following theorem, which has not been 
stated explicitly in the literature heretofore. 

Theorem 3.5 For the system x = Ax + Bu(t), 

if the homogeneous system is ASIL 
and has a Lyapunov function 
V = x^Qx with V = -x^Cx and 
1 1 u( t) | | <_ k 2 for all t, then the 
system is stable and the state is 
sure to enter a region defined by 

2k 2 IIQII ||B|| - -f- 

||x||<k, = (3.12) 

A 1 

where A, = A.„ (C) the minimum 

1 min 

eigenvalue of matrix C and e < 0 
with | £ | arbitrarily small 

Proof: The system is stable by definition 2.1 of Section I I C . Take 

as the system Lyapunov function V = x"*Qx then 
V = -x T Cx + 2 x T QBu 

^ -x T Cx + 2 1 1 x | | ||Q|| | | B 1 1 1 1 u 1 1 
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(3.12a) 



4 

also since 

x ] | | x | | 2 <_ x T Cx 

certainly <_ -X-j | |x| | 2 + 2k 2 | |x| | ||Q|| ||B|| 

To complete the proof we must show for all x such that | |x| | >_ k.| 

then V <_ e < 0 

When | |x | | >_ , then 

. . 2k HQII | | B | | - -f- 

IMI 1 L 

X 1 

and At I |x| I 2 > 2k 2 ||x|| ||Q|| | |B| | - ^i|iL 

Thus 0 > e X 1 1 |x| | 2 + 2k 2 | |x| | ||0|| ||B||>_tf 



It follows that 0 is negative when ||x|| >_ k-j , indicating that the 
systems response is leaving the region. 

The results of this theorem only give an upper bound for the region 
to which the state will converge in the steady state. Unfortunately , 
this bound depends upon the Lyapunov function chosen. However, the 
bound is determined without finding the solution to the system equations. 

Example 3 . Consider the stable system 



x 



A 



-1 

-1 

-1 

-1 





sin tot 



B 



1 

0 



x^Cx 



Relation (3.12a) is derived using the Reyleigh Quotient 
See Bellman [B3] p. 110. ||x 
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-1 ± j 



2 0 



MA) 




2 



and compute a steady state bound for the 
state of the system. 

Theorem 3.3, page 25 applies to the homogeneous dynamics, therefore a 
suitable Lyapunov function is 



V = x T x 

that is, Q = I 



V = 



C 



x T (A + A T )x 
1 6 ~ 

0 2 



-x T Cx 



and the homogeneous system is ASIL. The results of Theorem 3.5 are 
used to compute an upper bound for ||x|| in the steady state. The 
parameters for equation (3.12) are as follows taking for the matrix 

M = (m. .) ; ||M|| = [tr(MM T )]^ 2 (see Appendix A) 

* J 

1 1 u(t) 1 1 = 1 1 sin wt| | <_ k 2 = 1 

linn = i 1 1 b 1 1 = i 

x min< C > - 2 

e < 0 and | c | is arbitrarily small 



It follows that 



|x| | < k-|-' 



2kJ|Q| 



1 



It should be noted here that if C = X^I the computed bound k-| 
is independent of A, . This follows because -C = A T Q + QA. If for 
C = I the resulting Q is Q-j , then for C = A^I the resulting 0 is 
AiQ-|. Thus A-j cancels out in the expression for k^ . 
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The same theorem holds for discrete systems with appropriate 
changes in notation. 

It is important to note that, for linear systems, the forcing 
functions have no affect upon the asymptotic stability properties of 
the system. They do however, influence the behavior' of the system 
in the steady state. 



G. TRANSIENT ESTIMATION 

Since the level curves of V give, in some sense, a measure of 
distance from the origin, V may be used to estimate the convergence 
rate of asymptotically stable systems. The following is taken from 
Kalman and Bertram [Kl] with some minor additions. 

Consider the obvious inequality: 



¥(x,t) = V(x,t) < - n] V(x,t) (3.13) 

where n-j is the minimum of the ratio - in some region 

(radius r) of state space excluding the origin. In cases where the 
basic stability theorem applies (Theorem 3.1 page 21) we may take 5 



inf 

x 



Y,( 


X 


) 


e( 


X 


) 



0 < 



(3.14) 



Integrating equation (3. 13) observe that 



t 



/ 



t 



o 




< 



t 

/ n-j dt 



5 ^^(y), 0 < llyll < r> implies the greatest lower bound 

of f(y) for all y such that 0 < | |y| | <_ r. 
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- n l(t-t 0 ) 



(3.15) 



V(x(t),t) ^ V(x 0 ,t Q ) e 



The number, 

describing the convergence of the Lyapunov function following the 



Tjy , is an upper bound for the time constant 



system trajectory. If the system is linear the Lyapunov function is 

2 



a quadratic form and 



'1 



is the corresponding bound for the con- 



vergence of the system trajectory. 

Similarly we can write >_ -r^V where n 2 is a maximum of the 
ratio — y- , giving a lower bound for the system time constant. 

Kalman [Kl] has carried out the minimization to obtain n-j 
foi the linear time invariant case. These results are stated here. 

If the system x = Ax is ASIL and has a Lyapunov function V = x^Qx 
with 0 = -x^Cx we take 

'T 



'1 



ro in <x'Cx, V-l 



(3.16) 



The minimization leads to 

"i - W^' 1 ) < 3 - 17 > 

-Similarly 

"2 * W*' 1 * (3J8) 

Figure 1 illustrates the relationship of V(x(t)) and those functions 
formed by the n bounds. 

These convergence estimates depend upon the Lyapunov function 
chosen. However, there are certain cases where the estimates are 
identical to those obtained directly from the system matrix A. The 
following new result is stated as a theorem the proof of which uses 
a result indicated in Bellman [Bl] involving the eigenvalues of 
rational functions of several matrices, as follows: 
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V(x(t) ) 




Lyapunov Function and convergence bounds 
following the system trajectory x(t) 




FIGURE 2 

Effect of system parameters G on V. 
G minimizes V 
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Let A and B be commutative matrices and 
f(x,y) be a rational function. Then the 
characteristic roots of f(A,B) are 
) , where a; and b-j are the 
eristic roots of A and B in an 
ordered sequence which depends upon the 
matrices A and B. 

Theorem 3,6 If Theorem 3.3 of Section IIIC page 25 

anplies and Q is chosen as the identity 
matrix then -n-j = 2x - n (A) and the 
convergence estimate using Lyapunov 
functions is identical with that of the 
actual system. 

T 

The proof is immediate since if 0=1 then C = -(A + A ). 
Taking B = A^ in the above result and noting that a^(A) = x..(a"*"), 
we have that \.{A + A T ) = 2a i - (A). And the x ml * n ( C ) = 2x m .. n (A). 

The utility of this method is obviously not in estimating the 



f(a is bi 

charact 



convergence of linear time invarient systems since one eigenvalue prob- 
lem is replaced by another. Its utility in this thesis is in application 
to time varying systems as illustrated in the following example. 



Example 4. 



Let A(t) = 



-o4> 



o 

-< 2 4 > 



It is desired to estimate the convergence 
rate for the system described by the given A(t). 
Write this as 



A(t) = A - G(t) 




By applying the lemmas developed in Section IV C which follows, it can 
be shown that a suitable Lyapunov function is (see example 6 Section IV C) 



V = 



x T x 



Si = -x r (A + A T )x - x T (G(t) + G T (t))x (3.19) 
-Y = -x^Cx where -C = A + A**" 
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For this example 



Q = I 



C 



2 

0 



n l 




0 

4 




A . (Cl 
mi n 




2 



Theorem 3.6 page 35 applies also, hence 

n, = -2x . (A) = 2 

1 min' ' 



Suppose that \l is a function of a matrix of parameters G, as for 
example, the G(t) in equation (3.19). Now if for some G = G , 

V is a minimum with respect to G then the ratio — y— is also a 

minimum, hence, n-| , is maximum (see eq. (3.13)). It follows that if 
n-j is maximum for G = G then the system involving G converges 
faster than systems involving other G's. Figure (2) illustrates this 
concept of a minimum ratio -y- (G) which also appears in Section V 
in connection with the derivation of the optimal filter gain. In this 
case it is possible to choose G(t) so as to maximize the system con- 
vergence rate. 

Corresponding results to the foregoing are established by Kalman 
in [K2] for discrete systems. The main difference is that a recursive 
relation for V(x(n)) is established as 
V(x(n+1 ) ) < e' 2nT V(x(n)) 

where T is the sample period and n is the corresponding constant in 
the continuous case. Again the minimization leads to 



e -2"lT 



min 



-1 



) 
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H. CONCLUSION 



General results are difficult to obtain using the foregoing methods 
for time varying linear systems. Often the results give sufficient 
conditions only, and upper bounds which depend upon the chosen Lyapunov 
function. However, there is an advantage in that these results are 
always obtained easily without solving the system equation. 

The approach of restricting the Lyapunov function to be explicitly 
time invariant was considered and serves as an introduction for a unified 
approach to the special system problem considered in the next chapter. 
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IV. THE HOMOGENEOUS FILTER OBSERVER EQUATION 



A. INTRODUCTION 

The remainder of this thesis is concerned with the study of 
the basic recursive averaging filter-observer using the foregoing 
concepts of stability theory. 

The filter problem is usually formulated by considering a signal 
model driven by white noise and whose output measurements contain 
additive white noise. Such a model is described by 

x = Ax + w(t) (4.1 ) 

z = Hx + v ( t ) 

where A is an nxn state transition matrix, H is an mxn measurement 

matrix, w(t) and v(t) are n and m vectors respectively of uncorrelated 
white noise. 

Filtering is accomplished by recursively operating on the 
ovservations z to estimate the values of x. The estimates are 
denoted x. The filter dynamics are given by 

x = Ax + G(t)(z-Hx) (4.2) 

where G(t) is the nxm filter gain matrix the choice of which determines 
the filter characteristics. If the error, e = x-x, is formed and 
the expected value of the error squared (P = E[ee^]) is minimized 
to determine G(t) then (4.2) becomes the optimal integral squared 
error filter derived by Kalman and Bucy [K3]. If G(t) is determined 
by some other method, the filter is suboptimal (for example see S2). 

The optimal gain is given by 

G(t) = P(t)H T R _1 (4.3) 

where P(t) satisfies the matrix Ricatti equation 
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(4.4) 



P = AP(t) + P(t)A T -P(t)H T R _1 HP(t) + W 
and R - E[vv T ], U = E[ww T ]. 

Suppose we rewrite equation (4.2) so that it is in the form of 
equation (4.1). We obtain 

x = (A-G(t)H)x + G(t)z (4.5) 

This equation is interpreted as a system having a state transition 

matrix, A-G(t)H > with the input vector z. Thus, we are indeed 

taking the observations z and filtering them to provide the output 

A 

state vector x. 

Let us also find the equation for the error dynamics. Define 
the error, e = x-x, to be the difference between the actual states 
and the estimated values out of the filter. 

e = x-x = Ax+ w(t) - Ax - G(t)(z-Hx) (4.6) 

Noting from (4.1) that z = Hx + v(t) we write 

e = Ae - G(t)[Hx + v(t) - Hx] + w(t) (4.7) 

and e = (A-G(t)H)e - G(t)v(t) + w(t) (4.8) 

Observe that the homogeneous dynamics for the filter equation 
(4.5) and the error equation (4.8) are identical. Both are described 
by the state transition matrix A - G(t)H. Consequently the results 
of this section apply to both the output of' the filter and to the 
actual error. 

Now suppose that G(t) is a matrix of constants, G c . Then the 
filter transition matrix becomes A - G c H and the problem reduces 
to an ordinary linear time invariant case. Moreover, if G = I 
the filter-observer equation reduces to that form derived by 
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g 

Luenberger [LI] for observation of states. 

In discrete time problems the signal model and optimal filter 
equations becomes respectively: 



x(n+l ) 


= $(T)x(n) + w(n) 


(4.9) 


z(n) = 


Hx(n) + v(n) 




x(n) = 


4>(T)x(n-l) + G(n)[z(n) - H<i>(T)x(n-l )] 


(4.10) 


G(n) = 


P(n/n-l)H T [HP(n/n-l)H T + R]" 1 


(4.11) 


P(n/n) 


= P(n/n-l) - T G(n)HP(n/n-l ) 




P(n+l/n) 


= <t>P(n/n)i> + W 





where the notation P(n+l/n) is the estimated covariance of error 
at time n+1 given the previous n corrected estimates. 

For discrete systems, we find equations 

x(n+l ) = ( <t> ( T ) - G(n)H$(T) )x(n) + G(n)z(n) (4.12) 

e(n+l) = (4>(T) - G(n)H$(T) )e(n) - G(n)v(n) + w(n) (4.13) 
Again the unforced dynamics for the filter states and the estimation 
error are the same. 

Bona [B4] has considered the problem of choosing the eigenvalues 
of the matrix ( o (T ) - G(n)H$(T)) to give a desired convergence 
property (namely, exact estimate, of the observation in a finite 
number of sample periods for noiseless systems) leading to the 
specification of G(n). He has found that the observer can be made 
to converge as rapidly as desired. However, the observer is highly 
susceptible to error due to the presence of noise. 



Luenberger [LI] has shown by using linear transformation 
theory that unobserved output states of a free linear 
system may be observed by passing the observed output 
states through another linear system whose transition 
matrix is A-H, provided A-H is stable. 
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The filter gain matrix G(t) or G(n) has a twofold purpose in 



the filter observer equation. 

1. To provide observable outpus of unovserved states. 

2. To provide some programmed weighting of the actual 
observations in an attempt to compensate for 
errors introduced by noise. 

Further-properties -of G are studied~'in the following sections. 

B. FILTER STABILITY 

Stability of the optimal filter has been demonstrated in general 
by Kalman. A more general statement of the following theorem 
appears in Kalman and Bucy [K3] . 

Theorem 1 (Kalman) [ K3] . For a linear time invar- 

iant signal model also having stationary 
statistics, if the model is: 

(i.) completely observable 

(ii.) completely controlable 

then the optimal filter is asymptotically 

stable. 

It is evident that, even if the signal model equation (4.1) 
section IV. A is unstable, the filter may still be asymptotically 
stable. That is, the matrix (A-G(t)H) in equation (4.5) has 
its own stability properties. Thus, it is not necessary that the 
signal model be stable for the filter to be stable. 

Deyst and Price [Dl] have derived conditions for asymptotic 
stability of the discrete optimal filters by following the previous 
work of Kalman for the continuous optimal filter. Their method 
uses a time varying Lyapunov function namely choosing Q = P~"*(n) 
where P(n) is the discrete covariance of errors. Then»by applying 
the conditions for control abi 1 ity and observability, they found 
the desired bounding functions needed to apply the basic Lyapunov 
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stability theorem. The techniques developed in this paper, namely 
requiring the Lyapunov function to be time invariant, generally 
apply to a more limited class of signal models (i.e., asymptotically 
stable ones). However, the method is also much simpler and applies 
to suboptimal filters as well. Such a study leads to insight into 
what the filter gain must be to insure a stable filter and may in 
some cases be used to actually choose the G matrix requiring that 
the filter be stable. 

Generally, the technique proposed here may be outlined as 
follows : 

Assume that there is a Lyapunov function of the form V = e^Qe 
then compute 

V = e T (A T Q + QA)e - e T (H T G T (t)Q + QG(t)H)e (4.14) 

In some cases equation (4.14) may be used directly to determine 
the matrix Q. Whereas, in other cases, Q may be determined from 
-C = A T Q + QA by applying Theorem 3.2 page 22. If A is asymptot- 
ically stable then V of (4.14) becomes 

\l = -e T Ce - e T (H T G T (t)Q + QG(t)H)e (4.15) 

If equation (4.14) is negative definite for t> t then the 
filter is ASIL by theorem 3.1 on page 22. Negative definiteness of 
V may be shown for all t > t . Alternately it is possible to find a 
bounding function, y(||x||), as required by theorem 3.1. Obviously, 
if the second quadratic form in equation (4.15) is at least positive 
semi-definite for all t > t Q then, \l <_ -e^Ce and the conditions 
of theorem 3.1 are met. 
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The fact that the filter is ASIL means that considering the 
filter as a linear system free of input forcing functions, the 
error starting at any arbitrary initial value will tend to zero 
as t increases. 

This tells us nothing about the actual estimates or error. 
These questions must be considered with the input present which 
due to its stochastic nature, limits the asymptotic value of the 
effective error. Such consideration appears in Chapter V. 

For discrete systems, taking the same Lyapunov function as 



for the continuous case and using Theorem (3. 2D) page 24, we obtain: 





AV(e,n) = e T [$ T ( I - G(n)H) T Q(I - G(n)H)$ - Q]e 
AV(e,n) = -e T ce - e T ($ T QG(n)H$ + $ T H T G T (n)Q$)e 


(4.16) 




+e^$"^H^G^QGH4>e 


(4.17) 


and 


AV = -e T Ce - 2e T ($ T QG(n)H$)e + e T 4> T H T G T QGH$e 


(4.18) 


when 


-C = $ T Q$ - Q 





The following serve to illustrate the kinds of filter problems 
considered in this thesis. 

1 . Filter Gain Constant 

Case 1 We are given G(t) = G c and G c is a matrix 

of known constants. For this case we may 
take = e^Qe and rewrite equation (4.14) 
as 

V = e T ((A-G c H) T Q + Q(A-G c H))e (4.18a) 

Since A,H and G c are all known, we may 
apply Theorem (3.2) page 22 directly to 
determine Q. An algorithm for numerical 
generation of Q is derived in appendix B 
for the discrete case (Theorem 3. 2D). 
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Case 2 



Choice of G so that the filter is stable, 
c 

Now use equation (4.15) for V and determine 
Q based upon the signal model' by applying 
Theorem 3.2 page 22. This requires that the 
signal model be asymptotically stable if 
a suitable Q is to be found. 

The second quadratic form in equation (4.15) is now: 

e T (H T G T Q + QG H)e = e T De 

D = H T G r T Q + QG H (4.19) 

If G can be chosen so that the matrix D in (4.19) is 
c 

at least positive semi-definite(p.s .d. ) , the filter is sure 
to be stable since the first quadratic form in equation 
(4.15) is already negative definite. 



Example 5 



Consider again the system described in 
Example 3, page 30. Where 





-1 1 




1 


A = 




B = 






-1 -1 




0 



Q = I 



C = 



i — 

2 

0 



0 

2 



Let H = [1 0] 



In case 2, it is desired to choose G c so 
that the matrix, D, in (4.19) is at least 
positive semi-definite, which is a suffi- 
cient condition to insure that the filter 
be stable. In order that H and G„ be 
conformable, G c must be of the form 
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then 



G H = 
c 



91 0 

9 2 0 

and from (4.19) 

^9] q 

D = 1 92 

g 2 0 

The condition that D be p.s.d. is that its 

principal minors be non-negative [ref. B2], 

That is 

2g-| 1 0 



2 

J 2 - 



-go > o 



Obviously g2 = 0 and g.| >_ 0 satisfy the 
conditions. Therefore, any filter with 



0 



and gi >_ 0 



will be ASIL. 

Now suppose we are given 

r 

0 

Find the convergence bound o-| 
V = e T (M T q + QM)e 



From (4.18a) 



where 



M = A-G H = 
c 



-2 1 
-1 -1 



We may compute the Lyapunov function Q from 
theorem 3.2 for any p.d. matrix C such that 



-C = M Q + QM 
But in this case 



/|T - an( j we ma y use 
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theorem 3.3, page 25. Thus f 

j +4 0 

Q = I C = -(M + M ) = 

0 +2 

The convergence bound, r^, is given by (3.17) 

"1 - W")" 1 ) x min^ ’ 2 

2. Special Time Varying Gain 

Let the filter gain matrix be time varying such that 

lim G(t) = 0 (4.20) 

t^°° 

the null matrix. Here we have a problem very similar to that 
discussed in section III E, and those comments apply here. Also 
if the matrix 

(H T G T (t)Q + QG(t)H) (4.21) 

in equation (4.15) is at least positive semi-definite for all 
t > t Q , then we may not need the requirement of (4.20). Sufficient 
conditions for the matrix (4.21) to be positive semi-definite are 
developed in section IV C which follows. 

3. The Optimal Filter 

Substituting the optimal gain G(t) = P(t)H^R”^ into 

equation (4.15) we obtain 

V = e T (A T Q + QA)e - e T (H T R _1 HP(t)Q + QP(t)H T R _1 H)e (4.22) 
If the signal model is asymptotically stable, Q may be determined 
for any p.d. matrix C as before and (4.22) becomes 

V = -e T (C + H T R' 1 HP(t)Q + QP(t)H T R _1 H)e (4.23) 

Now the condition for stability from (4.23) is that 

f(y,t) = y T (C + H T R _1 HP(t)Q + QP(t)H T K _1 H)y > 0 (4.24) 

for all t. Moreover, if 



+4 0 

0 +2 
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min f(y,t) > 0 
t 

then the condition will hold for all t. From ordinary calculus 
a minimum occurs for t = t-j if 

df( *'V = 0 and , o 

~~ 

Taking the derivative of equation (4.24) yields 

f(y,t) = y T (H T R _1 Hf>(t)Q + Qf'(t)H T R" 1 H)y = 0 (4.25) 

One way for (4.25) to be zero is for 

P(t) = 0 (4.26) 

Kalman [K3] has shown that condition (4.26) yields an asymptotic 
or steady state value for P(t) in the limit as t approaches infinity 
by solving the Ricatti equation (4.4) page 39 for P with f 1 = 0. 
However, for this condition, lim P(t) is also zero, hence, the 
foregoing does not in general lead to the desired minimum. Only 
in the trivial case where P(t) is diagonal does this method lead 
to the desired stability condition, since the diagonal elements of 
P(t) are monotonic decreasing to the limit. 

The elements of the matrix in (4.25) are, in general, linear 
combinations of all the elements of P(t). There may be some 
t = t-| which will make (4.25) zero other than in the limit as 
t approaches infinity. However, ^ is generally difficult to 
find, even for a specific example. Often it is necessary to compute 
the trajectory for all elements of P(t); then, using this, study 
the properties of equation (4.24) at various values of time in an 
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effort to determine t^ . Besides being difficult, the method ceases 
to be useful, since, once P(t) is known for all t, complete knowledge 
of the filter is at hand. 

Further research is required in the area of finding a 
suitable bounding function, y(||e||), for equations (4.22 or 4.23). 

For suboptimal filters one might choose G(t) in equation 
(4.15) in order to establish the stability condition. The next 
section attempts to build some insight into the nature of G(t) 
required to insure a stable filter. 

C. STABILITY CONSTRAINTS ON G(t) 

In the last section we have found for the filter error dynamics 
that if V = e^Qe is a Lyapunov function then 

\l = -e T Ce - e T (H T G T (t)Q + QG(t)H)e (4.27) 

V = -e T Ce - 2e T (QG(t)h)e (4.27a) 

where -C = A^Q + QA and Q may be determined if the signal 

model x = Ax is ASIL. Equation (4.27a) is written as shown 

for convenience even though the second term is not a true quadratic 
form ( i . e . , the weighting matrix is not necessarily symmetric). 

In this section some sufficient conditions are found which insure 
that V in (4.27) be negative for all t. If these sufficient conditions 
hold for a given G(t) then the filter is shown to be ASIL. On the 
other hand, one may choose G(t) in order to make the conditions apply, 
thus designing a filter which is known to be ASIL. 

In applying the basic stability theorem of section III B, we 
must find a positive function y such that 
V(e,t) <_ - y( | | x | | ) for all t 
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Equivalently, we may use a quadratic form for the bounding 
function of tf, that is, we must find a positive definite matrix, 

D, such that: 

\7(e,t) <_ - e T De = -y (4.28) 

Applying this to equation (4.27) we require: 

e T [C + H T G T (t)Q + QG(t)H]e > e T De (4.29) 

At this point we introduce inequality notation for symmetric 
matrices. Namely, A > B means that (A-B) is positive semi -definite. 
If equality is removed, then (A-B) is positive definite. Writing 
relation (4.29) in this notation: 

C + H T G T (t)Q + QG(t)H ^ D for all t (4.30) 

The matrix D is important if one is to use the technique 
described in section 3G to estimate transient response of the filter. 
The form y= e^De replaces in the ratio 7/V used to determine 
n-| hence estimating the convergence rate of the filter. Large D 
means that the gradient of ^ is steep, hence, the convergence bound 
is better since the actual filter converges faster than the estimate 
given by n -| . Thus, if G(t) is known, then it is desired to choose 
the largest D which will satisfy (4.30). On the other hand, if we 
wish to design the filter to converge faster than some given bound, 
then D is given and (4.30) provides a constraining relation for 
G(t). 

If we rewrite (4.30) and apply the properties of norm to both 
sides, we get another bounding relation for G(t). 

H T G T (t)Q + QG(t)H > D-C (4.31) 



49 



| I H T g‘ ( t)Q| I + | |QG(t)H| | >_ D-C 

2IIQII I |G(t)| | ||H||> IID-CII 



1 G(t) 1 1 > 1 




D-C 


| for all t. (4.32) 


2 


IQI 


nwr 



Intuitively equality in (4.31) and (4.32) should result in a 
minimum condition of the left hand sides for all t. That is in 



(4.32) 



min | | G(t) | | = 


D-C 


1 M 


t 2 


IQI 


l|H|| 



(4.32a) 



A similar interpretation of the equality statement for (4.31) 
is not so easily conceived. 

The basic requirement for stability given by (4.30) or (4.31) 
is that 



min x . n (H T G T (t)Q + QG(t)H - D + C) = 0 (4.33) 

t 

Application of (4.33) is very difficult. However sufficient conditions 
to insure that (4.33) does hold, may be found by applying an important 



theorem involving bounds on the eigenvalues of a matrix. The Hadamard- 
Gerschgorin theorem is stated here as found in Bodewiq [B3] . 



Theorem (Hadamard-Gerschqorin) 



The eigenvalues of the complex matrix 
A = [a. .] lie inside the closed domain 

* J 

G consisting of all circles 
with centers a. and radii r. where 

- • S - J 



K.j(i = 1 .,n) 



If we apply the H-G theorem to a real symmetric matrix and 
require that a.-... _> / |a..| for all i then we can get a sufficient 



condition test for positive semi-definiteness. 
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Lemma 1 . Let A be a real symmetric matrix and if 

a ii - ^I|a.j . | for all i, then A is positive 
semi-definite. 

The proof is straight forward application of the H-G theorem. 
Since A is real symmetric its eigenvalues (a(A)) must be real. 

Then applying the H-G theorem, the eigenvalues of A must lie between 
min ( a -jj _ r -j ) an d max + ) where r 4 = / ,|a 4 _. | as 






illustrated in Figure 3a. But the hypothesis requires 
implying 



a ii l r i 



a^ - r^ >_ 0 for all i 

In particular, min (a-. - r.) >_ 0, 
i 11 ' 1 

which implies that the eigenvalues of A must be non-negative. Hence, 
A must be positive semi-definite and the lemma is proved. 

Two lemmas will now be proven which will provide a method by 
which we may determine the best bounding matrix D. 

Lemma 2 . Let: 

a. M(t) be a real symmetric matrix of elements 

m i j (t > 

b. m 1i (t) - 

c. a.j(t) >_ 0 for all i and t 
Then M(t) is positive semi-definite for all t. 

The proof is immediate from Lemma 1. Figure 3b illustrates 
the a . (t) . 
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Hadamard-Gerschgorin Theorem applied to 
a real symmetric matrix 



A = 




r ii 




Im[x] 
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As a result of lemma 2, if from equation (4.31) we set 



M(t) = H T G T (t)Q + QG(t)H 
then M(t) > D-C for all t 



(4.34) 

(4.35) 



If inf a . (t) = 0 for any i then then condition for stability 

t 1 

is obtained by choosing D = C. However, if a.(t) > 0 for all t, 
then we should be able to choose some D > C giving rise to a tighter 
restriction on G(t). 

Lemma 3 . In Lemma 2 let: 

a. a. (t) = min a- (t) for all t 

. i 

b. a^(t 1 ) = min a^(t) > 0 



The proof of Lemma 3 is by construction. Choose t-j as specified 
in conditions (a) and (b), thus insuring that the smallest eigen- 
value of M(t) will be positive for all t. For the matrix (M(t) - M(t-j)), 
use (b) in Lemma 2 to compute condition (c) in Lemma 3, and the 
lemma is proved. 

It is apparent that the conditions in Lemma 3 may be used to 
insure that equation (4.33) holds. That is, by choosing M(t) as 
in (4.34) and the matrix M(t-j) as in Lemma 3, we may rewrite 
equation (4.33) as an inequal ity since Lemma 3 only specifies a 
boundary for the eigenvalues of M(t). The result is 



c. 




Then M(t) >_M(t^) with equality only at t = t^ . 



if 



min Nn-inW*) - D + O) >_ 0 




(4.36) 

(4.36a) 



then 



53 



Now if Lenina 3 applies, then we are assured that (4.36a) is true 
for some choice of the matrix D. However, we may also apply Lemma 3 
directly to (4.31). This results in 

M(t) > D - C (4.37) 

Since C is given and D may be computed directly, by replacing the 
left side of (4.37) with its minimum, M(t-|) , to establish the equality 
M(t 1 ) = D - C (4.38) 

This follows, since by Lenina 3, M(t) >_M(t^). 

On the other hand, if D is given, as established by some con- 
vergence bound, n-| , then the conditions in Lemma 3 may be applied to 
(4.37), thus, establishing constraining relations on the elements of 
M(t). 

We are really interested in the matrix G(t). From (4.34) it is 
obvious that every element of M(t) will in general be a linear combi- 
nation of all the elements of G(t). This follows, since H and Q 
are determined by the signal model and are constant matrices. Thus, 
in general 



M(t) = Tm kl ] = Cb,, g,. <t)) kl ] 



(4.39) 



ij 



This may not be difficult to apply to specific examples since H and 
Q are often very simple matrices and many of the b. . = 0. 

The following illustration refers to Example 4, page 35. 

Example 6 . As in Example 4, consider the system described by 

the transition matrix 



-d + 4-) 



A(t) = 



-(2 + 4 -) 
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This is equivalent to a filtering problem 
where the signal model is described as 
follows: 



The system: % - 



-1 0 
0 -2 



x + 



w 



w. 



and observation 



z = 



where A 



-1 0 
0 -2 



1 


0 


X 


+ 


■»r 


0 


1 






V 2 



H = I 



and w-j , w 2 , v-j , v^ are white Gaussian noise 
inputs. The homogeneous filter dynamics are 
given in (4.8) as 



e = (A - HG(t) ) e 
where 



A(t) e 



G(t) = 



1 



t 

0 



0 

J_ 

t 



By using Lemma 2 we will derive the Lyapunov 
function given in equation (3.19) of Example 4. 

First, assume a Lyapunov function of the form 

V = e T Qe 

then by equation (4.15) 

V = -e T Ce - e T (H T G T (t)Q + QG(t)H) e 

T n 

^ -e De 

where Q is determined from -C = A^"o + OA. But in this case AA^ = a"*A 



and we may use Theorem 3.3, page 25 to find Q and C . The result is 
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Q = I 



C 



-(A T + A) 



2 

0 



0 

4 



Now in order to apply Lemma 2 we find from equation (4.34) 
M(t) = G T (t) + G(t) = 



This follows since Q = H = I . 

Then by Lemma 2 we must have 

2 

“l = a 2 = ~t~ — 0 for all t. 

Obviously, a ^ 0 for any t >_ t Q > 0. Since ^ a(t) = 0 in the 
limit as t approaches infinity, the condition for stability is ob- 
tained by choosing D = C . That is, -e^Ce is a suitable bounding 
function for V . Hence, \( is negative definite and the filter is 
ASIL by Theorem 3.1, page 21. 

Lemma 3 applies to problems where the o^(t) in Lemma 2 are greater 
than zero for all t. Consider the slight extension of the above exam- 

0 



pie by letting 



G(t) = 



t 

0 



+ 1 




By equation (4.34) 



M(t) = 



2 + 



0 



And in Lemma 3 



_ , . 2 
“i ■ 2 + ~r 



4 + 



“2 ’ 4 + T 
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Take t-j in the limit as t approaches infinity. Then from conditions 
(a) and (b) in Lemma 3: 
a{t,) = 2 > 0 

It follows that 

2 0 

0 4 



Mft,) = 



and condition (c) Lemma 3 reduces trivially to Lemma 2. 
(4.38) we may compute the matrix D 
M^) = D - C 

"4 0 



D = M(t 1 ) + C 
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Using equation 



By using the results of Section III G we may compute a convergence 
bound for the filter employing the given G(t). 

The problem of using the foregoing Lemmas to obtain constraints on 
the elements of G(t) giving rise to a class of suboptimal filters is 
deferred to Chapter VI. 

The foregoing Lemmas apply to discrete systems simply by replacing 
t with n . Now M(n) is specified from equation (2D) in Section 4B. 
Namely: M(n) = 4* T QG(n)H$ + 4>W(n)Q<f> 



- $ T H T G T (n)QG(n)H$ 



(4.40) 



However, (4.40) may be simplified. This follows since 4> is non- 
singular and appears as a congruent transformation on the inclosed 
matrix and does not affect the sign definiteness of the inclosed matrix 
[ref. B2]. So in effect we may take 



M(n) = QG(n)H + H T G T (n)Q - H T G T (n)QG(n)H 



(4.41) 
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D. CONCLUSIONS 



In this section we studied the homogeneous dynamics of the basic 
filter-observer equation (4.5) and its error dynamics equation (4.8). 
Since these homogeneous dynamics are identical, the results of this 
section apply to both. 

The basic approach in this section is to assume a time invariant 
Lyapunov function and compute its time derivative, 0 , subsequently 
showing that \l is negative definite for a given filter gain matrix 
G(t). This was found difficult to do in general terms. Often it 
is possible to determine the quadratic form for V based solely 
upon the signal model. This requires that the model be ASIL. 

f 

It is desirable to choose the Lyapunov function based only upon the 
signal model if one wishes to compare the convergence of two 
different filters designed for the same signal model by using the 
estimation techniques described in section III G. Under these 
conditions the Lyapunov function for the two filters is the same, 
but the associated V is modified according to the particular filter 
(see equation (4.15)) . 

The techniques developed in this chapter cannot be readily 
applied to the optimal filter. This is left as an area for addi- 
tional research. 

Several constraining relations between the bounding y function 
required for stability and characteristics of the gain matrix were 
developed in section IV C. In general these relations are difficult 
to apply, but in some simple cases they lead to a method of deter- 
mining G(t) to insure a stable filter which converges faster than 
some desired bound. 
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V. THE FORCED FILTER EQUATION 



A. INTRODUCTION 

In this chapter the forcing functions are included in the 
filter dynamics. From the study of linear forced systems in 
section III F, recall that the forcing functions do not affect 
the asymptotic stability of the system dynamics. They do, however, 
determine the steady state properties of the system. Such is also 
the case for the filter-observer system, and these properties 
are investigated in this chapter. 

The error dynamics for the filter given in equation (4.8) 
are repeated here 

e = (A-G(t)H) e - G(t) v(t) + w(t) (5.1) 

As in chapter IV, assume a Lyapunov function of the form 

V = e T Qe (5.2) 

then by direct computation using (5.1) 

• «T T • 

v = e Qe + e Qe 

= e T (A T Q + QA) e - e T (H T G T (t)Q + QG(t)H) e 
+ 2w T (t)Qe - 2v T (t)G T (t)Qe 

V = -e T Ce - e T (H T G T (t)Q + QG(t)H) e 

+ 2w T (t)Qe - 2v T (t)G T (t)Qe (5.3) 

where Q is determined, as before, for any given p.d. matrix C by 
the signal model which is required to be ASIL. 

It should be noted that the only difference between equation 
(5.3) and equation (4.15) is the addition of two terms involving 
the forcing functions w(t) and v(t). However, this particular 
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application is different because the forcing functions are random 
vectors. To consider this aspect we use the expectation operator^ 
and rewrite equation (5.3). Since Q is completely determined for 
any choice of positive definite C, we write this equation in terms 
of Q. Applying the matrix identity (5.4) to (5.3) 

x T My = tr[Myx T ] = tr[yx T M] (5.4) 

results in 

V = tr[Q(A-G(t)H)ee T + ee T (A-G(t)H) T Q 

+ 2Qew T (t) - 2ev T (t)G T (t)Q] (5.5) 

In order to remove the random nature from this problem we take 
the expectation of V and V. The forcing functions w(t) and v(t) 
are uncorrelated white Gaussian random processes. Consequently, 
the estimate error is also a Gaussian random process. We define 
the error covariance matrix P = E[ee^]. Thus the Lyapunov 

function (5.2) becomes: 

E[V] = trE[Qee T ] = tr(QP) (5.6) 

The expectation operation is done over an ensemble of Lyapunov 
functions each of which is the result of one experiment. Taking 
the expectation and applying (5.4) to (5.5) results in 
E[tf] = tr[Q(A-G(t)H)P + P(A-G(t)H) T Q] 

+ 2QE[ew T (t)] - 2E[ev T (t)]G T (t)Q] (5.7) 



^The expectation operator is defined as E[x] = / xf(x) dx 
where f(x) is the density function of x 
[ref. Papoulis, PI] 
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To finish this problem we must first determine what E[ew ] 

and E[ev^] are. Equation (5.1) has a solution of the form 

t t 

e(t) = <t> ( t-, t )e(t ) + / $(t,x)w(x)dx - / ^ ( t , x )G(x ) v (t ) dx 

L o L o (5.8) 

But w(t), v(t) and e(t Q ) are uncorrelated so that 

E[w T w] = E[v T v] = E[w T e(t 0 )] - E[v T e(t 0 )] - 0 (5.9) 



Also define 



T 



E[w(t)w (t)] = W6 (t-T) 

EO(Ov T (t)] = RS (E-t) 
where a(t-t) = \ £ 

Now using (5.8), form ew T . 

e(t)w^(t) = $(t,t )e(t )w T (t) + / $(t,x)w(x)w^(t)dx 

00 t 

o 



(5.10) 



- / 4>(t,x)G(x)v(x)w T (t)dx 

t 



(5.11) 



Taking the expectation of (5.11), interchanging the order of 
integration and expectation and using (5. 9), (5. 10) 

E[e(t)w T (t)] = 0 + / t *(t,x)Wfi(t-x)dx - 0 



= 1/2 [W] (5.12) 

g 

This follows from the sifting property of the Delta function . 
Similarly we can find: 



E[e(t)v'(t)] = -1/2 [G(t)R] 



(5.13) 



o 

The Delta function sifting property applied to an integration 
limit at the location of the Delta function is used here, 
b 

(i.e., / f(x)6(x-b)dx = 1/2 [f (b)] 

a 
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Substituting (5.12) and (5.13) into (5.7) yields 

E[V] = tr[Q(A-G(t)H)P + P(A-G(t)H) T Q + QW + G(t)RG T (t)Q] 

(5.14) 

Note that if Q = I in (5.14) the result becomes the error 
dynamics developed by others. Specifically Athens and Tsi [Al] 
arrived at (5.14) in their derivation of the optimal filter and 
Sims and Melsa [S2] also arrived at this result for specific optimal 
filters . 

For discrete systems a difference equation for E[V(n)] analogous 
to equation (5.14) is easily obtained by using equation (4.13) 
to form 

E[e T Qe] = tr[QP(n)] (5.15) 

similarly as for the continuous case. The result for P(n) is 

P(n) = (I-G(n)H)P'(n)(I-G(n)H) T + G(n)RG T (n) (5.16) 

P 1 (n) = 4>P(n-l )$ T + W 

where R is the covariance matrix for the measurement noise and W 
is for the noise driving the signal model. A derivation of equation 
(5.16) also appears in Sorenson [S3]. 

B. DERIVATION OF THE OPTIMUM G(t) 

In effect the stochastic Lyapunov function (5.6) for the filter- 
observer is a linear combination of all the covariance of errors 
between the estimates (outputs of filter) and the actual signal 
being estimated. The covariance of error gives a measure of how 
well the filter should be expected to perform. As such, the Lyapunov 
function (5.6) gives a scalar measure of the performance of the 
filter. It is desirable that the covariance of error be minimum 
for all t. Thus the Lyapunov function (5.6) may be used as a cost 
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criterion to be minimized. The optimum filter is then formed by 
choosing G(t) such that E[v] is minimized for all t. If (5.14) 
is negative for all t >_ t , then 

E[V] < 0 (5.17) 

and the filter is ASIL since E[v] > 0 for all t. Moreover, if 
E[V] is negative, then E[V] is assured to be a minimum if the 
magnitude of its time derivative is maximum for all t. The optimum 
G(t) can be derived as follows, using the concept of a gradient 
matrix. The following formulas for the partial derivative of a 
scalar function (tr) of several matrices with respect to one of 
the matrices can be proven (Athans and Tsi [A1 ] ) . 

8 tr(AB T ) = 8 tr ( BA T ) = A (5.18) 

8 B 8B 

8 tr(ACA T ) = 2AC (5.19) 

3 A 

where C is symmetric. 

Taking the partial derivitive of (5.14) with respect to G(t) 
using (5.18) and (5.19) yields 

8 E[V] = 0 = - 2PH T Q + 2G(t)RQ 

8G 

and solving for 

G(t) = PH T R _1 (5.20) 

This is the identical result derived by Kalman and others 
[A1 , K3, S3] for the optimal filter. To complete the derivation 
of the optimal filter we must show that (5.20) results in a minimum 
for (5.14). Note the second partial derivative of (5.14) with 
respect to G, by applying (5.18) and the fact that R and Q are 
symmetric. 
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(5.21) 



If QR = RQ 



(5.22) 



R and Q commute, then the right hand side of (5.21) is positive 
definite since Q and R are positive definite. This follows from 
the result indicated in Bellman [Bl] given on page 33 of chapter III. 
Now if (5.21) is positive definite then (5.20) does result in a 
minimum for (5.14). But if (5.17) holds, then the magnitude of 
E [ty] is maximum (for G in (5.20)), hence E[V] in (5.6) is a minimum. 

Moreover, the ratio 



is also maximum, and the following important conclusion is reached. 

(l) The filter using the G(t) given 
by (5.20), is assured of having 
the most rapid convergence rate 
to the minimum covariance of error. 

This follows from the discussion on page 36 of chapter III, section G. 

The foregoing derivation strictly holds only if conditions 
(5.17) and (5.22) are true. However, the result indicated in 
(1) above does give engineering insight into the nature of optimal 
filters and to desired properties of suboptimal filters. That is, 
it is possible to make an engineering trade-off between convergence 
rate and steady state covariance of error in the design of suboptimal 
filters which are much easier to implement than the Kalman filter. 

C. CONCLUSION 

The concept of a Lyapunov function for random variables was 
introduced for the forced filter dynamics. A derivation of the 



-E[V] 

E[V] 



(5.23) 
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optimal filter gain matrix was given leading to the conclusion that 
the optimal filter is the most rapid converging to the minimum 
covariance of error. As a consequence, 5 it is possible in the design 
of suboptimal filters to make an engineering trade-off -between 
convergence rate and steady state estimation error. 
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VI. AN EXAMPLE OF A SUBOPTIMAL FILTER 



A. INTRODUCTION 



Consider the discrete form of a system described in phase variable 
form with only one state being observed. Namely: 



x = Ax 



y = Hx 



where 



A = 



0 1 0 



0 



0 1 
-a 



H = [1 0 ... 0] 

( 6 . 1 ) 



It is assumed that A is stable and H is a 1 x n matrix, and 
that $(T) has been found as described in Section I I A where T is 
the sampling period. Note the transition matrix for the discrete 
filter (4.12) 

x(n+l) = ($(T) - G(n) H $(T))x(n) + G(n)z(n) 



If H is 1 x n, then G(n) must be n x 1 in order for the indicated 
multiplication to be conformable. This is indeed fortunate, since G(n) 
is chosen to insure that the filter is stable by application of the 
lemmas in Section IV C, and these lemmas only give n constraining 
relations for the elements of G. 

Proceeding, we determine the system Lyapunov function V = e"*Qe 
for aV = -e^Ce by finding Q such that 't^Qo - 0 = -C for a given 
positive definite C. This may be easily done using the algorithm 
developed in Appendix B. 

Now, in order to insure stability, the matrix M(n) in (4.41) must 
be positive semi -definite in accordance with Lemmas 1 and 2 developed in 
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Section IV C. The argument, (n) .will be dropped for convenience 
in the following and the notation used for Q is: 



G ~ ^ij^ n x n 



and for G = 



n x 1 



Repeating equation (4.41) here 

M = QGH + H T G T Q - H T G T QGH 



( 6 . 2 ) 



Multiplying (6.2) out, noting that: 
G = [g 1 g 2 . . . g n ] T 

GH = [G i f 0] n x n 



H T G T = 



r 

0 



QGH = 



f b ij g d 



n 

E b . q . 
•| nj “j 







x n 



G QG is in fact a quadratic form in the vector G and 



H T H = 



1 0 . . . 0 

0 0 



then H T G T QGH = ( z b. . g. g.) H T H 

i j 1 J 1 J 

Substituting the results, (6.3) through (6.4), into (6.2) 

( l Vj - ^ l VV-) 1 V.i • • • I 



z b 2j gj 



(6.3) 



(6.4) 



M 



(6.5) 



Now the 



of Lemma 2 are: 



= 0 - 



n 

l b ii 9i 

j=l 1J J 



i = 2,3, ... ,n 



According to Lemma 2, a . > 0 however, it is evident that for 



l - 



i ^2 we must choose ou = 0 since they can never be positive. This 
results in n - 1 linear equations in n unknowns depicted by: 



n 

z 

j-1 



b. . g . = 0 

ij a J 



i = 2,3, ...,n 



( 6 . 6 ) 



Choosing to solve these n - 1 equations for g^ ... g n in terms of 
g.| , rewrite equations (6.6) in matrix form 



1— — 








— 


b 22 b 23 * * • b 2n 
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" b 21 


b 32 




9 3 




' b 31 


b n2 ‘ ’ b nn 




9n 




b nl 









1 


- 



9 1 (6.7) 



The square matrix in (6.7) must have an inverse, since it is a 

O 

principal minor of Q which is positive definite. Thus, (6.7) may in 
fact be solved for g 2 ... g p in terms of g-|. The remaining relation 
necessary to determine g^ is given by of Lemma 2. 



a l 



1 b ljSj 



b. . g . g . - 0 

,• U 1 J 



j=l u i»J 

But must be non-negative, giving: 
n 

E b,. g . > E b g g 

j=l IJ 3 i , j 1J 1 J 



( 6 . 8 ) 



(6.9) 



g 

The determinants of the principal minors of a positive definite 
matrix are positive [ref. B2]. 
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Using equation (6,7) and relation (6.9) a family of G matrices may 
be determined. Since some Is zero, the stability condition as 
discussed in Section IV C following Lemma 2 is obtained by choosing 
D = C. Referring to Section III G the convergence bound for this filter 
is: 

e ' 2ni ■ W 1 *' 1 * ■ W*' 1 ) t6 - 10) 

But, (6.10) is the convergence bound for the signal model since C 
was chosen for the signal model. This means that the transient behavior 
of the filter is assured to be no worse than that of the signal model. 

Unfortunately, the foregoing Lyapunov technique developed cannot be 
used directly for optimal filters involving the phase variable form of 
signal model. In particular, aV , equation (4.18), has not been shown 
to be negative definite for the optimal filter in this example. How- 
ever, a reasonably good comparison of the suboptimal filter design and 
the optimal filter may be made by computing tr(QP(n)) , equation (5.15) 
for both filters using the Q determined from the signal model. This 
follows because the Lyapunov function (5.15) is a scalar function of all 
the covariance of error for the filter. 

B. NUMERICAL EXAMPLE 

In order to summarize the foregoing results, a simple third order 
numerical example is given. 

A suboptimal fi lter-observer is to be designed for the given signal 
model. The filter is required to be ASIL. 

The signal model is described as follows. Given the signal dynam- 
ics described by the continuous time matrix in phase variable form 
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A = 



0 1 0 
0 0 1 
-15 -17 -12 



( 6 . 11 ) 



the discrete version of this signal model for a sample period of T = .15 
is driven by a white Gaussian random sequence whose covariance matrix is 
given as 



W = 



0 

0 

0 



0 0 

0 0 

0- 0.0086 



( 6 . 12 ) 



The corresponding discrete system for (6.11) was found for T = .15 
by using a digital computer. For this, a common NPS subroutine called 



AT 



.994 


.143 


.0065 


.0977 


.884 


.0653 


.979 


-1.207 


.101 



"PHIDEL" was used which essentially computes the series for e 
resulted in: 



*(.15) = 



Only one state of the signal model is observed and 
H = [1 0 0] 

The observations are contaminated by additive white noise whose 
covariance matrix is 



This 



(6.13) 



(6.14) 



R = 



.5 

0 

0 



0 

.5 

0 



0 

0 



(6.15) 



This completes the description of the signal model. 

The first step in the design of a suboptimal filter is to find a 
Lyapunov function for the signal model. 



70 



Assuming that C = I the Q matrix of the Lyapunov function for 
this system was found using the algorithm developed in Appendix B. For 
this problem the algorithm converged in 42 iterations to a stopping 
criterion of 



max 

ij 



d. .(k) 
ij 



.0005 



where d. .(k) is an element of the k 
Appendix B for further explanation. 

The resulting Q is: 



th 



Q = 



13.9 

7.56 

0.357 



7.56 

13.6 

0.789 



0.357 

0.789 

1.08 



correction matrix. 



(6.16) 

See 



(6.17) 



In order to design the filter, G must have the n x 1 form 



G = 



91 

92 

9 3 



(6.18) . 



and equation (6.7) must be solved. 

(6.7) and solving for gg and g^ 

r ~ir 



13.6 0.789 



92 



0.789 1.08 




Substituting from (6.17) into 
results in 



-7.56 

-0.357 



(6.19) 



1 

CM 

CP 




0.768 -0.0561 




-7.56 
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-0.561 0.968 




-0.357 






— > __ 




_ 



9l 



9 2 = -0.56 g-j 
g^ = 0.055 g.| 



( 6 . 20 ) 
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The constraining relation (6.9) is used to obtain a restriction 
on g-j. For this example (6.9) becomes 

* 2 

13.9 g.| + 7.56 g 2 + 0.357 g 3 >_ 13.9 g^ + 2(7.56) g-^ 

+ 2(0.357) g.jg 3 
+ 13.6 g 2 2 
+ 2(0.789) g 2 g 3 
+ 1.08 g 3 2 

Substituting the relations (6.20) and reducing this equation yields 
9.69 g-j >_ 9.69 g-j 2 

g 1 <_ 1 for g 1 > 0 (6.21) 

Thus if g^ is constrained by equation (6.21) for all n , and 

g 2 ,g 3 are determined by (6.20), the filter is sure to be stable by 

its construction using the sufficient conditions of the Lyapunov method 

in conjunction with the Hadamard-Gerschgorin Theorem. 

The simplest kind of filter to implement is one with constant gains. 

A simple time varying gain which satisfies (6.21) is: 

g-. — , for n — 1,2,3,... 

1 n 

As discussed in Section VI A, a comparison of these filters with 
the optimal filter can be made by plotting tr(QP(n)) for each filter. 

The covariance of error matrix P(n) was computed on a digital 
computer for the optimal filter and for the following two suboptimal 
gains established by (6.20) and (6.21). 



~1 


i 


~l/n 


-.56 


G 2 = 


-.56/n 


.055 




.055/n 


(_ _J 




— M 
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Computer subroutines for these calculations appear in the Computer 
Program Section. The subroutine "GAIN" was used for the optimal filter 
and "SUBCOV" for the filter gains in (6.22). Then using the Lyapunov 
matrix, Q , in (6.17) 

E[V] = tr(QP(n)) (6.23) 

was plotted on the same axes for all three filters as shown in Figure 4. 

It is interesting to note that all three filters are stable, and 
the optimal filter does converge faster than the suboptimal filters as 
predicted in Chapter V, Section B. The optimal filter has its greatest 
advantage in the first one or two sample periods. That' is, the correc- 
ted covariance of error (P) is changing more rapidly here than for the 
suboptimal filters. 

C. CONCLUSIONS 

The main developments in this thesis show that it is feasible to 
design suboptimal filters which do perform nearly as well as the Kalman 
filter. This is of great advantage because suboptimal filters often 
lead to a much simpler implementation. Such filters are designed so 
that they will be ASIL and converge faster than a known bound. The 
design is accomplished by using Lyapunov functions and requires that the 
signal model be ASIL. 

There are two areas brought out by this thesis which require further 
research. 

1. Find a suitable bounding function for the time derivative of 
the Lyapunov function (4.22 or 4.23) for the optimal filters. 

2. Find bounding relations for the covariance of error matrix, P , 
similar to those developed in Theorem 3.5 for the system state. 
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The two most important properties of the filter-observer system are 
its convergence rate and its ability to discriminate against noise. By 
combining any results obtained for 1. and 2. above, with the techniques 
developed in Chapter IV of this thesis, a design procedure for sub- 
optimal filters could be formulated which would insure that the filter 
be stable and also give information about its two most important prop- 
erties. Hence, such designs could be formulated and compared with 
other such designs without the necessity of simulating each filter. 
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(20 units per inch) 




n (10 units per inch) 



FIGURE 4 

Lyapunov Function for Various 
Filter-Observers 



Curve 


Filter gain 




OPT 


OPTIMAL 




SUB1 


G ] = (1 -.56 


.055) T 


SUB2 


G 2 = (1/n -.56/n 


.055/n) 
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APPENDIX A 
NORMS 



A norm is a scalar function of several variables which satisfy 
the following conditions. The norm is denoted ||*|| . 

i) | |X| | > 0 and | |X| | = 0 iff X = 0 

ii) ||wX|| = |oj| ||X|| where |u| is the magnitude 

of a complex scalar 

iii) | |X + Y| | _< | | X | | + | | Y| | Minkowski inequality 



Possible norms 


for an n vector x 




INI 


n 

= ^ l*il 

1-1 1 


ii*ii - i“l„ 1*, i 


IMI 


- ( s lx/ ) ,/2 - 


(tr(xx T )) ,/2 


11*11 


- U|x/) ,/P 


8 

V 

D- 

V 


Possible norms 


for a general matrix 


A 


MAN 


= [tr(AA T )] 1/2 


l|A|| = e | a. . | 

i > J J 


IIAII 


- T " i a i -s i 

1 j=i 1J 


l|A|| - J | a . . | 2 
i j J 


MAN 


= sup 1 1 Ax 1 1 
1 1*11 = 1 





Other properties — Schwarz inequality 

||AB|| < ||A|| x ||B|| ||AX|| < ||A|| x ||x|| 



Further details may be found in Ref. K1 , SI, T1 . 
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APPENDIX B 

Lyapunov Function for Linear Discrete System 

Davison and Man [D3] have derived an algorithm to numerically 
generate a Lyapunov function for linear time invariant systems 
given a negative definite quadratic form for \l . Here a similar 
algorithm is developed for'use with linear discrete systems. 

In section III C it is shown that for the discrete system: 



x(n+l) = A q x ( n ) (1) 

a positive definite Lyapunov function of the form: 

V = x T Qx (2) 

exists with: 

AV = -x T Cx (3) 

being negative definite. If the system (1) is stable^ Q is 
determined for any positive definite € from the relation 

-C = A d QA d - Q (4) 



In general it is difficult to solve equation (4) for Q even 
given that C = I let alone some arbitrary positive definite C. 

So another approach is taken and assumes that the system described 
by (1 ) is ASIL. 

Proceeding, we write the Lyapunov function following the 

system trajectory described by (1) as: 

V ( n ) = V(0) + V(l) - V(0) + V( 2) - V(l) + ... + V(n) - V(n-l) 

since AV(k) = V(k) - V(k-l) 

n 

then V(n) = V(0) + s AV(i) (5) 

i=0 
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Now let n approach infinity then V(n) approaches zero by 
assumption that (1) is ASIL. Then (5) becomes: 

V(0) = ■ 

Using equation (3) 

V(0) = 

Noting that the solution to (1) is: x q and that 

v(0) = xjqx 0 

We rewrite equation (6): 

x o <> x o - j x o < A D> TcA D x o 



Z AV(i) 
i=0 



°° T 

Z x i C Xi 

i=0 1 1 



( 6 ) 



which implies that: 

Q = z (Ah T CAj 

i =0 u u 



this summation will converge if equation (1) is stable, since 

| A ( A d ) | < 1 and will not converge if (1) is unstable. 
The proposed algorithm as stated here is easily implemented 
on any computer system having appropriate matrix subroutines. 



D k " A D D k - 1 A D 



^k+1 



= Qix + d l 



k = 0,1, 2, 3,... 



with initial conditions: 



% • c 

and termination criterion: 




The matrix, D, is the correction matrix and a simple conver- 
gence criterion is to terminate the iteration when the magnitude 
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of all elements of the correction matrix are less than some desired 



amount. That is: 

max | d . . | <_ e where D = (d . .) 
ij J 

Subroutine "LYAP", given in the computer program section following, 
was used to generate the Lyapunov Q matrix for C = I in the third 
order example given in chapter VI. The result is given in equation 
( 6 . 17 ). 
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c 

c 



THIS SUBROtJ T I NE COMPUTES THE OPTIMUM GAIN MATRIX AND 
THE ERROR COVARIANCE 



SUBROUTINE GAI N(PKK,PKKM1 ,Q,R,PHI,H,N,M,G,HI, ND,MD, LO) 
ODIMENSI ON PKK ( 12,12) ,0(12,12) ,H(12 , 12) 

1 ,G( 12,12) ,R( 12,12 ) ,H I (12,12) ,HT( 12.12 ) ,TEMP( 12, 12), 

2, TEMP 1( 12.12) , PHI ( 12,12) ,PHIT( 12, 12 ) , PKKMl ( 12, 12 ) 
3TEMP2( 12,12) 

C G(K) = P(K/K-1)*HT*(H*P(K/K-1 )*HT + R)**(-l) 

CALL TRANS(H,N,N,HT,ND,MD) 

CALL PRCO( PKKMl ,HT,N,N,N,TEMP,ND,MO,LD) 

CALL RR OD(H, TEMP, N,N,N, TEMPI, ND, MO, LO) 

CALL A 00 (TEMPI , R ,N,N , TEMPI , NO , MO) 

CALL REC!P(M, 0.000001 , TEMPI ,TEMP2, KER, MO) 

IF (KEP-2) 101,110,101 

110 WPITE(6,111) 

111 FORMAT (5HKFR=2) 

101 CALL PR00(TEMP,TEMP2 ,N,N,N,G,ND,MD,LD) 

C NOTE HERE PKK(I,J) = P(K/K) WHERE 

C P(K/K) = (I-G(K)*H)*P(K/K-1) 

CALL PPOD(G,H,N,N,N,TEMP,NO,MD,LD) 

DO 108 1=1, N 
00 108 J= 1 , N 

103 TEMP( I,J)=-TEMP(I ,J) 

CALL ADD (HI , TEMP, N,N, TEMP, NO, MO) 

CALL PROD (TEMP, PKKMl , N , N ,N , PKK, ND, MD, L 0 ) 

C NO T E HERE PKKMl (I, J) = P(K/K-1) WHERE 

C P(K/K-1)= PHI*P(K-1/K-1)*PHIT + Q 

CALL TRANS(PHI , N , N , PH I T , NO, MO) 

CALL PROD (PKK, PHI T ,N , N , N,TEMP , NO, MO, LO ) 

CALL PROD(PHI , TEMP, N,N,N, TEMPI ,ND,MO,LO) 

CALL AOD( TEMPI, 0,N,N, PKKMl ,ND, MO) 

RETURN 

END 



SUBROUTINE SUBC0V(PKM1 , W . R , PHI , H , N, M, G , H I , ND, MO, PK ) 
DIMENSION PK(12 ,12) ,PKM1 (12,12) ,W ( 12 , 12 ) , R ( 12 , 1 2 ) , 
CH( 12,12 ),G( 12.12 ),HI( 12,12 ) ,TEM (12, 12), 0(12, 12), 
#PHI( 12.1Z),PHIT(12,12) 

PK=( I-GH l PKMi ( I-GHJT+GRGT 

CALL PROD(G,H,N,N,N,TEM,ND,ND,NO) 

CALL SUB(HI ,TEM,N,N,TEM,ND,NDr 
CALL TRANS(TEM,N,N,PHIT,ND,ND) 

CALL PP00(TEM,PKM1 , N ,N, N, 0 , NO, NO, ND ) 

CALL PPOD(D,PHIT,N,N,N,PK,ND,NO,ND) 

CALL TRAN$(G.N,N,D,NO,ND) 

CALL PROD(G,R,N,N,N,TEM,ND,NO,ND) 

CALL PROD(TEM,D,N,N,N,PHIT , NO, NO, NO) 

CALL AOD(PK,PHlT,N,N,PK , N, ND) 

PKM1=PHI*PK*PHIT+W 

CALL TRANS(PHI , N , N , PHI T , ND , NO ) 

CALL PROD ( PK , PHI T, N,N,N, TEM,ND,NO, ND) 

CALL PROD (PH I ,tEM,N,N,N,PKM1 ,ND,ND,ND) 

CALL ADD(PKM1,W,N,N,PKMI,ND,ND) 

RETURN 

END 
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SUBROUTINE LYAP( sphi » OL » N » NO , MP , C ) 

DIMENSION SPHI (12,12) ,QL ( 1 2 , 1 2 ) ,C ( 1 2 , 1 2 ) , PH IT ( 1 2 , 12 ) , 
«D1( 12,12) ,02(12,12) ,X2( 900) 

C LYAPUNOV FUNCTION 
OP 2004 1=1, N 
DO 2004 J=1,N 
0 2(I,J)=C(I,JJ 

2004 QL( I, J)=C( I , J) 

CALL TRANS(SPHI , N , N , PH I T ,MD , NO ) 

DP 2005 K=1 , 500 

CALL PPOD( D H I T ,D2,N,N,N,01 ,NO,ND,NO) 

CALL PROP (01 ,$PHI ,N,N,N,02 ,NO,Nn,ND) 

CALL ADO ( OL ,02 » N» N , QL »N0 ,ND) 

TST=0. 

00 2006 1=1, N 
00 20 C 6 J= 1 , N 
TT p =A B c ( 02 ( I ,J) ) 

2006 CONTI NUF 

IF( TTP.GT.TST) TST=TIP 
IF( Tcj.LE. 0.0005) GO TO 2007 

2005 CONTINUE 

2007 WPITE (6, 2008) K 

2008 FOP MA T ( * ITTERATION K=',I3) 

RETURN 

END 
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derivation shows that the optimal filter converges most rapidly. A design of 
a suboptimal filter for one class of signal models is given based solely upon 
stability constraints. 
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